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Abstract 

A method for classifying types of brain activity in 
magnetoencephalographic (MEG) signals is proposed. 
Sources of abnormal cortical activity are localized by 
performing a generalized spectral analysis in the space of 
Fourier coefficients of the expansions of recorded signals in 
adaptive orthogonal bases. The basic principles of the 
method are discussed, and the results of its application to 
actual MEG records are presented for functional brain 
mapping in normal and pathological states. 
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Introduction 

Magnetoencephalography (MEG) is a recently 
developed experimental technique for studying 
human neural activity, mapping functional brain areas, 
and diagnosing various disorders. Research in this 
area is motivated by the possibility of noninvasive 
imaging of cortical and subcortical processes. 

Until recently, the technique was difficult to apply 
because of the extreme weakness of both spontaneous 
and evoked neuromagnetic fields generated by current 
sources. Accordingly, the measuring apparatus to be 
employed must meet particularly stringent 
requirements. 

Despite technical difficulties and the high cost of the 
required apparatus, this technique can be used as 
complementary or alternative to electroencephalo- 
graphy because magnetic field is much less perturbed 
by intracranial inhomogeneities and integumentary 
tissues than electric field, and the accuracy of source 


localization can therefore be substantially improved 
without detailed knowledge of the intracranial 
structure. 

Data to Study 

The input experimental data used in the present study 
were obtained by means of a 148-channel Magnes 2500 
WH neuromagnetometer (a highly sensitive SQUID- 
based on magnetic field measurement system) at the 
New York University School of Medicine. The system 
was housed in a magnetically shielded room. 

The analysis presented in this paper is focused on 
recognition and classification of patterns of magnetic 
field distribution over the subject's scalp for localizing 
neuromagnetic activity sources. This problem arises 
from the complex variability of an MEG signal during 
the measurement time (approximately 10-20 min), 
which includes alternating intervals of abnormal and 
normal activity. Intervals of abnormal and normal 
activity can easily be identified visually when they 
differ by the mean amplitude of the signal. Otherwise, 
visual inspection cannot be reliably used to identify 
intervals corresponding to different activity modes, 
and the dominant type of activity can be determined 
by analyzing the spatial field pattern, which strongly 
correlates with the number, location, and orientation 
of the current sources responsible for magnetic field 
generation. 

The proposed approach was validated against data 
derived from a reference group of healthy individuals 
and a group of patients suffering from Parkinson's 
syndrome. The signal was a 148-dimensional time 
series measured at 148 points on the head surface at a 
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sampling rate of 500 Hz per channel. 


pass, high-pass and band-pass filtering. 


The data analysis can be divided into the following 
tasks: 


■ detection of a signal associated with a 
particular type of brain activity (e.g., response 
to auditory, visual, tactile, or other stimulation; 
generation of Parkinson's tremor or auditory 
hallucinations; etc.); 

■ selection of instants for sampling the spatial 
field distribution over the head surface to be 
used as input in solving the inverse problem of 
current source localization; 

■ solution of the inverse problem of source 
localization in the cases of abnormal and 
normal activity; 

■ the use of fMRI data as physiological 
constraints. 

Measurements with Phantom Sources 

For the setting of algorithms for inverse problem 
solving the preliminary measurements were carried 
out. Here the data used have been obtained in 
measurements with phantom spherical conductor and 
were generated by 10 fixed dipoles of variable 
moments. Fig. 1 demonstrates this case. On the right 
figure points correspond to sensors positions, black 
and white are for positive and negative flux of the 
magnetic induction respectively. 



FIGl. TEST DIPOLES 1 - 5 AND THEIR ANTIPODES (LEFT) AND 
EXPERIMENTAL DATA SET FOR DIPOLE 1 (RIGHT) 


An advantage of such problem formulation was in a 
priori knowledge of dipole positions. Thus, we can 
compare results of computer calculations in inverse 
problem with real dipole coordinates. 

MEG Data Pre-processing 

The reprocessing refers to the reading, segmenting 
around interesting events such as triggers, temporal 
filtering and optionally re-referencing of the data set. 
In the work preprocessing implements a Fourier 
transform filter, a Butterworth and a FIR filter for low- 


Detection of artefacts can be done byusing manual or 
automatic artefact detection and visual identifying 
independent components (IC A). There are two ways to 
remove artefacts: 1) rejecting trial or the piece of trial 
containing the artefact; 2) subtracting the spatio- 
temporal contribution of the artefact from the 
data.ICA is most widely used for cleaning from 
cardiac, eye movements and blinks artefacts. 

The application of transcranial magnetic stimulation 
(TMS) pulses during the data acquisition poses some 
specific challenges. This is an occurring of TMS 
artefact right after TMS pulse (pulseartefact, 
ringing/step response artefact, cranial muscle artefact) 
and recharging artefact in500ms after stimulation 
onset. There is a general trend for increase in artefact 
magnitude with TMS intensity. The following artefacts 
removal methods were used: removal large muscle 
artefacts from TMS by ICA (Korhonen et al., 2011), 
adapted recording system with software artefact 
correction (Litvak et al., 2007), off-line removal of 
TMS-induced artefacts on human magnetoencephalo- 
graphy by Kalman filter, based on cross-correlation 
coefficients of ICA (LeeU. et al., 2013). Thentime 
average ERP/ERF for all trials is computed 
andbaseline correction can be applied. 


Frequency Analysis 

It can be shown by means of spectral routines that 
oscillatory components contained in the MEG signal 
power changes regarding the system status. For this 
purpose frequency or time-frequency analysis based 
on Fourier bases and wavelets were used. 


As an example, monaural audible stimulus was 
applied. Pulsing frequency was 7 Hz, while the 
response was localized at frequencies near 10 and 20 
Hz. Fig. 2 shows this situation on the wavelet diagram 
(Makhortykh 2 ), 2009). 



FIG.2. HAAR WAVELET DIAGRAM OF MEG RESPONSE FOR 7 
HZ PULSING AUDIBLE STIMULATION 
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The Direct Problem 

Consider the dipole current density J generated in an 
active region G of conductive brain tissue 
characterized by conductivity a. Assuming a constant 
permeability p = po, we write the quasi-static Maxwell 
equations as 

E - -V (p , 

V B = 0 , (1) 

J = J p + crE , 

where <p is the electric potential and J is the total 
current density. The magnetic field is calculated by the 
Biot- Savart law: 

B(r) = ^JJ(r)xAzllj v '. (2 ) 

4 n | r — r | 

The variable of integration Y f is dipole radius- vector. 

Here we assume that the conductive region is a sphere 
and the source is a current dipole j(r) = Q£(r-r 0 ) with 
dipole moment Q and position vector of a source r 0 . 

According to (Sarvas, 1987), the magnetic induction on 
a homogeneous sphere enclosing the dipole can be 
calculated as 

B(r)=^F 2 (FQxr 0 -Qxr 0 -rV J P) (3) 

4 n 

It is clear that the outer magnetic field generated by a 
radially oriented source vanishes. 

Therefore, the inverse problem is ill posed: thesource 
characteristics can be determined only up to an 
arbitrary radial component of the current dipole 
moment. Note also that the contribution of volume 
currents to the magnetic induction given by Eq. (3) 
does not depend on the conductivity a. 

The Inverse Problem 

The locations and orientations of the sources of 
spontaneous activity are determined by using the 
magnetic field values at the measurement points to 
solve the inverse problem for Eq. (3). For the adopted 
model of brain activity and head geometry, the inverse 
problem can generally be reformulated as a nonlinear 
optimization problem for calculation of the locations 
and moments of the set of dipoles that provides the 
best approximation of the measured MEG data with 
respect to RMS deviation. The inverse problem is 
difficult solve primarily because it is ill posed: there 
exist infinitely many solutions consistent with the 
experimental data. The computational complexity of 


the nonlinear optimization problem is explained by 
the existence of multiple local extrema of the objective 
function, particularly when the number of dipole 
sources is large. 

There also exist technical difficulties due to the 
weakness of both spontaneous and evoked magnetic 
fields generated by brain current sources against a 
high background noise level in MEG measurements. 

To classify the types of neuromagnetic activity and 
identify abnormal intervals in MEG records, we use 
the feature space spanned by the spectral expansion 
coefficients of the total magnetic field on a sphere, 
which are invariant under rotation. The natural basis 
for describing the spatial field distribution is the set of 
spherical harmonics. Relevant features are extracted at 
the training stage of the algorithm by estimating the 
variability of the expansion coefficients for the class of 
objects under analysis. The upper bound for the 
dimension of the feature space depends on the length 
N of the Fourier expansion. There exists a simple 
analytical relationship between the expansion 
coefficients of functions whose arguments are related 
by an SO (2) transformation (Vilenkin, 1968), which 
makes it possible to develop a fast search procedure 
for functions in a particular class. 

Source Analysis 

Source estimation comprises two steps: 1) estimation 
of the potential or field distribution for a known 
source and for a known model of the head is referred 
to as direct modelling; 2) estimation of the unknown 
sources corresponding to the measured MEG is 
referred to as inverse modelling. 

The head model specifies how currents that are 
generated by sources in the brain, e.g. dipoles, are 
influenced by the tissue properties and how these 
result in externally measureable MEG fields. For MEG 
the following head models are used: analytical single 
sphere model; local spheres model for MEG, one 
sphere per channel; realistically shaped single shell 
approximation, based on the implementation from 
Guido Nolte; magnetic dipole in an infinite vacuum. 
In the earliest studies, head models with simple 
geometries and homogeneous parameters were used, 
permitting many simplifications in the computation. 
When using models with spherical symmetry, solving 
the direct problem can be reduced to evaluating an 
analytic expression in Eq. (3). 

To reconstruct the location and the time-course or 
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spectral content of a source in the brain, various 
source-localization methods are available. Here 
calculations with audible stimulation data and 
Parkinson patients' records used MRIAN software 
(Ustinin, Makhortykh et al., 2002). It is rather efficient 
for one- and two- dipoles problem with usage of 
structural restrictions by means of fMRI tomograms. 

Results of spontaneous and induced sources 
localizations are presented by Fig. 3. It shows the 
response on the monaural audible stimulus in 
auditory cortex: 10 Hz source is in the contra lateral as 
long as 20 Hz source (see Fig. 2) is in the ipsilateral 
hemispheres. Fig. 3 b) demonstrates one of the 
localized source of pathological activity for the case of 
parkinsonian patient. Full description is brought out 
in (Makhortykh 1 ), 2009). 



FIG.3. SOURCES LOCALIZATION (CIRCLES WITH RODS) FOR: A) 
INDUCED ACTIVITY (AUDIBLE STIMULUS RESPONSE) AND B) 
SPONTANEOUS PARKINSONIAN SIGNAL. 


Conclusions 

Analysis of MEG records and localization of brain 
biomagnetic activity are the bases of the 
cerebralprocesses functional mapping. Data 
preprocessing by means of generalized spectral 
method (Dedus, Makhortykh et al.,) allows to 
removing noise, highlighting the desired signal and in 
many cases determines the success of the solution. It 
uses a set of orthogonal bases so that spectral 
description can be scaled down to fit a particular class 
of signals. In many cases correct source localization 


and model parameters estimates can be obtained only 
after projection of the signal on the proper basic 
function (as on Fig. 2). Also expansion coefficients can 
be used as features for the system control and signal 
type recognition. 
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